clear all;
close all;
clc
H = 19100e3;
Rz = 6400e3;
Tobr = 12*60*60; 
Tmod = Tobr/2;
wobr = 2*pi/Tobr;
c_light = 3e8;
f_0 = 1602e6;

dt = 5*60;
t = 0:dt:Tmod;
rr = [0; Rz];
  rsv = (H+Rz)*...
      [cos(wobr*t);
      sin(wobr*t)];
  rv = rsv -repmat(rr,1,length(t));
  x_rv=rv(1, :);
  y_rv=rv(2, :);
  D = sqrt(rv(1, :).^2 + rv(2, :).^2);
  ind = find(rv(2, :) > 0 );

V = diff (D)/dt;
df = -V/c_light*f_0;
figure;
subplot(2,1,1)
plot(t(2:end), df, t(ind), df(ind - 1));
xlabel('t,sec')
ylabel('\Delta f, Hz')
subplot(2,1,2)
hist(df(ind-1), 10)
xlabel('Delta f, Hz')
 

